new_immigration_bo <- left_join(wall_u_attempt, immigration_data, by = c("ddyad", "year"))
View(new_immigration_bo)
write.csv(new_immigration_bo, "new_immigration_data.csv")
View(grav_vars)
new_immigration_bo2 <- left_join(new_immigration_bo, grav_vars, by = c("ddyad", "year"))
write.csv(new_immigration_bo2, "new_immigration_data2.csv")
write.csv(new_immigration_bo2, "new_immigration_data2.csv")
setwd("/Users/bailee/Documents/Replication Files for Weathering the Storm")
setwd("/Users/bailee/Documents/Rep)Files_WTS")
setwd("~/Documents/Replication Files for Weathering the Storm")
source(ReplicationFile_Conjoint_JASA.R)
library(MASS)
library(reshape2)
library(reshape)
library(countrycode)
library(states)
library(pltesim)
library(readstata13)
library("plyr")
library("foreign")
library("cjoint")
install.packages("FindIt")
library(FindIt)
library(MASS)
library(reshape2)
library(reshape)
library(countrycode)
library(states)
library(pltesim)
library(readstata13)
library(plyr)
library(foreign)
library(cjoint)
install.packages(FindIt)
library(FindIt)
setwd("~/Documents/Rep_Files_WTS")
conjoint<- read.csv("conjoint.csv")
conjoint$Gov <- factor(conjoint$Gov, ordered = TRUE,
levels=c("Democracy", "Anocracy", "Autocracy"))
conjoint$Pow <- factor(conjoint$Power, ordered = TRUE,
levels =c("top 15 most powerful", "top 25 most powerful", "top 50 most powerful"))
conjoint$Reg<- factor(conjoint$Reg, ordered = TRUE,
levels = c("Middle East North Africa", "subSaharan Africa",
"Asia", "South America", "Europe"))
conjoint$All <- factor(conjoint4All, ordered = TRUE,
levels = c("Defensive Alliance", "Offensive Alliance", "Neutrality Pact"))
conjoint$Rep <- factor(conjoint$Rep, ordered = TRUE,
levels = c("Consistently fulfills alliance commitments",
"Inconsistently fulfills alliance commitments",
"Rarely fulfills alliance commitments"))
conjoint$HR <- factor(conjoint$HR, ordered = TRUE,
levels = c("Excellent human rights record",
"Mediocre human rights record",
"Routinely violates human rights"))
conjoint$select <- factor(conjoint$selected, ordered = FALSE,
levels = c("0", "1"))
conjoint$Gov <- factor(conjoint$Gov, ordered = TRUE,
levels=c("Democracy", "Anocracy", "Autocracy"))
conjoint$Pow <- factor(conjoint$Power, ordered = TRUE,
levels =c("top 15 most powerful", "top 25 most powerful", "top 50 most powerful"))
conjoint$Reg<- factor(conjoint$Reg, ordered = TRUE,
levels = c("Middle East North Africa", "subSaharan Africa",
"Asia", "South America", "Europe"))
conjoint$All <- factor(conjoint$All, ordered = TRUE,
levels = c("Defensive Alliance", "Offensive Alliance", "Neutrality Pact"))
conjoint$Rep <- factor(conjoint$Rep, ordered = TRUE,
levels = c("Consistently fulfills alliance commitments",
"Inconsistently fulfills alliance commitments",
"Rarely fulfills alliance commitments"))
conjoint$HR <- factor(conjoint$HR, ordered = TRUE,
levels = c("Excellent human rights record",
"Mediocre human rights record",
"Routinely violates human rights"))
conjoint$select <- factor(conjoint$selected, ordered = FALSE,
levels = c("0", "1"))
View(conjoint)
conjoint$Gov <- factor(conjoint$Gov, ordered = TRUE,
levels=c("Democracy", "Anocracy", "Autocracy"))
conjoint$Pow <- factor(conjoint$Pow, ordered = TRUE,
levels =c("top 15 most powerful", "top 25 most powerful", "top 50 most powerful"))
conjoint$Reg<- factor(conjoint$Reg, ordered = TRUE,
levels = c("Middle East North Africa", "subSaharan Africa",
"Asia", "South America", "Europe"))
conjoint$All <- factor(conjoint$All, ordered = TRUE,
levels = c("Defensive Alliance", "Offensive Alliance", "Neutrality Pact"))
conjoint$Rep <- factor(conjoint$Rep, ordered = TRUE,
levels = c("Consistently fulfills alliance commitments",
"Inconsistently fulfills alliance commitments",
"Rarely fulfills alliance commitments"))
conjoint$HR <- factor(conjoint$HR, ordered = TRUE,
levels = c("Excellent human rights record",
"Mediocre human rights record",
"Routinely violates human rights"))
conjoint$select <- factor(conjoint$selected, ordered = FALSE,
levels = c("0", "1"))
fit10 <- CausalANOVA(formula = selected ~ Rep + Pow +
Reg + All + HR + Gov,
data = conjoint, diff=TRUE, pair.id = conjoint$case_no,
cluster = conjoint$random_id, nway=2, select.prob = TRUE,
seed = 45869, boot = 5000)
library(FindIt)
fit10 <- CausalANOVA(formula = selected ~ Rep + Pow +
Reg + All + HR + Gov,
data = conjoint, diff=TRUE, pair.id = conjoint$case_no,
cluster = conjoint$random_id, nway=2, select.prob = TRUE,
seed = 45869, boot = 5000)
fit10 <- CausalANOVA(formula = select ~ Rep + Pow +
Reg + All + HR + Gov,
data = conjoint, diff=TRUE, pair.id = conjoint$case_no,
cluster = conjoint$random_id, nway=2, select.prob = TRUE,
seed = 45869, boot = 5000)
fit4 <- FindIt(model.treat= selected ~ Gov + Pow +
Reg + All + HR + Rep,
nway=3, treat.type="multiple", data = conjoint)
print("This model has finally finished!!!! Yay!!!")
full.FindIt <- function(object){
## Make the full factorial design matrix as the target population.
comb.list <- apply(object$treat.orig,2,unique)
full <- as.data.frame(expand.grid(comb.list))
invisible(full)
}
full <- full.FindIt(fit4) #Produces full factorial design matrix
pred.dat <- predict(fit4, newdata=full, sort=FALSE)$data #Predicted potential outcomes
print("Predicted potential outcomes complete! Hurray!")
#Save model
save(pred.dat, file="Resolve Conjoint SVM Pred Dat.RData")
INT <- function (object, target.data, column, dist = "target", base,
sort = TRUE, compare = FALSE, order = 2)
{
if (all(class(object) != "FindIt")) {
warning("the class of object needs to be FindIt.")
}
Restriction <- NULL
measure <- "TIE"
if (compare == TRUE) {
if (missing(order)) {
warning("Need to specify the order")
}
}
coefs <- object$coefs.orig
names(coefs) <- c("Intercept", object$name.t)
outcome.type <- object$type
if (dist == "unique") {
data <- cbind(object$y.orig, object$treat.orig)
colnames(data)[-1] <- colnames(object$treat.orig)
data <- unique(data, MARGIN = 1)
}
if (dist == "sample") {
data <- cbind(object$y.orig, object$treat.orig)
colnames(data)[-1] <- colnames(object$treat.orig)
}
if (dist == "target") {
data <- target.data
print("Using Data representing the target population.")
}
else {
print("Using the Data used to fit the model.")
}
if (!missing(column) & !missing(base)) {
for (i in 1:length(column)) {
ind.column <- which(colnames(data) == column[i])
data[, ind.column] <- relevel(data[, ind.column],
base[i])
}
}
if (!missing(column)) {
ind <- c()
for (i in 1:length(column)) {
ind[i] <- which(colnames(data) == column[i])
}
column <- column[order(ind, decreasing = FALSE)]
}
base <- 1
Marginal <- function(data, base, Restriction) {
Effect.list <- list()
Name.list <- list()
range.mar <- c()
for (i in 1:(ncol(data) - 1)) {
columnM <- colnames(data)[i + 1]
A <- tapply(data[, 1], data[, columnM], mean, simplify = FALSE)
Effect1 <- unlist(A)
if (base == "min") {
base <- which(Effect1 == min(Effect1))
}
TEffect <- Effect1 - Effect1[base]
range.mar[i] <- max(TEffect) - min(TEffect)
Effect.list[[i]] <- TEffect
Name.list[[i]] <- paste(colnames(data)[(i + 1)],
names(TEffect), sep = "_")
}
names(range.mar) <- colnames(data)[2:ncol(data)]
name <- unlist(Name.list)
Treatment.Effect1 <- unlist(Effect.list)
Treatment.Effect1 <- as.data.frame(Treatment.Effect1)
colnames(Treatment.Effect1) <- "AMTE"
rownames(Treatment.Effect1) <- name
output <- list(Treatment.Effect1 = Treatment.Effect1,
range = range.mar)
return(output)
}
if (compare == FALSE & missing(column)) {
Treatment.Effect1 <- Marginal(data = data, base = base,
Restriction = Restriction)$Treatment.Effect1
Treatment.Effect <- Treatment.Effect1
}
if (compare == TRUE & order == 1) {
Marginal.range <- Marginal(data = data, base = base,
Restriction = Restriction)$range
print("Range of Marginal Effects")
return(Marginal.range)
}
comb.prem <- function(data, column, Restriction, order) {
A <- tapply(data[, 1], data[, column], mean, simplify = FALSE)
A2 <- tapply(data[, 1], data[, column], mean, simplify = TRUE)
A[is.na(A2)] <- NA
Effect1 <- unlist(A)
if (base == "min") {
base <- which(Effect1 == min(Effect1))
}
Comb.Effect <- Effect1 - Effect1[base]
Treatment.Effect <- cbind(Comb.Effect, expand.grid(dimnames(A)))
Treatment.Effect <- na.omit(Treatment.Effect)
Comb.Effect <- Treatment.Effect[, 1]
Combination.name <- Treatment.Effect[, -1]
Treatment.Effect.print <- Treatment.Effect
for (j in 2:ncol(Treatment.Effect)) {
Treatment.Effect[, j] <- paste(colnames(Treatment.Effect)[j],
Treatment.Effect[, j], sep = "_")
}
Treatment.Effect1 <- Marginal(data = data, base = base,
Restriction = Restriction)$Treatment.Effect1
Main.comb <- c()
for (i in 1:nrow(Treatment.Effect)) {
main <- c()
for (j in 1:(ncol(Treatment.Effect) - 1)) {
main[j] <- Treatment.Effect1[rownames(Treatment.Effect1) ==
Treatment.Effect[i, (j + 1)], 1]
}
Main.comb[i] <- sum(main)
}
Sum.Mar.Effect <- Main.comb - Main.comb[base]
if (order == 2) {
TIE <- round(Comb.Effect - Sum.Mar.Effect, digits = 8)
}
if (order == 3) {
data.column.three <- data[, column]
Two.way.comb.list <- list()
for (j in 1:3) {
column.com.three <- c()
Two.way.comb <- list()
column.com.three <- colnames(data.column.three)[c(combn(length(data.column.three),
2)[, j])]
A.Three <- tapply(data[, 1], data[, column.com.three],
mean, simplify = FALSE)
A2.Three <- tapply(data[, 1], data[, column.com.three],
mean, simplify = TRUE)
A.Three[is.na(A2.Three)] <- NA
EffectTwo <- unlist(A.Three)
Two.way.comb <- EffectTwo - EffectTwo[base]
Two.way.comb2 <- cbind(Two.way.comb, expand.grid(dimnames(A.Three)))
for (k in 2:3) {
Two.way.comb2[, k] <- paste(colnames(Two.way.comb2)[k],
Two.way.comb2[, k], sep = "_")
}
colnames(Two.way.comb2) <- c("Two.way.comb",
"first", "second")
Two.way.comb.list[[j]] <- na.omit(Two.way.comb2)
}
Two.way.comb.Final <- do.call(rbind, Two.way.comb.list)
Two.Way.sum <- c()
for (i in 1:nrow(Treatment.Effect)) {
two.sum <- c()
for (j in 1:nrow(Two.way.comb.Final)) {
ind.two <- sum(is.element(Two.way.comb.Final[j,
2:3], Treatment.Effect[i, 2:4]))
if (ind.two == 2) {
two.sum[j] <- Two.way.comb.Final[j, 1]
}
else {
two.sum[j] <- 0
}
}
Two.Way.sum[i] <- sum(two.sum)
}
Sum.Two.Effect <- Two.Way.sum - Two.Way.sum[base]
TIE <- round(Comb.Effect - Sum.Two.Effect + Sum.Mar.Effect,
digits = 8)
}
range.change <- max(TIE) - min(TIE)
Order.data <- cbind(Sum.Mar.Effect, Comb.Effect)
Order.data <- Order.data[order(Order.data[, 1], decreasing = TRUE),
]
Order.data <- as.data.frame(Order.data)
Order.data$order.main <- seq(1:nrow(Order.data))
Order.comb <- Order.data[order(Order.data[, 2], decreasing = TRUE),
"order.main"]
Order.diff <- Order.comb - seq(1:nrow(Order.data))
max.order.diff <- max(abs(Order.diff))
return(list(A = A, Treatment.Effect = Treatment.Effect,
Comb.Effect = Comb.Effect, Combination.name = Combination.name,
Sum.Mar.Effect = Sum.Mar.Effect, TIE = TIE, range.change = range.change,
max.order.diff = max.order.diff))
}
if (compare == FALSE) {
if (!missing(column)) {
X <- comb.prem(data = data, column = column, Restriction = Restriction,
order = order)
A <- X$A
Treatment.Effect <- X$Treatment.Effect
Comb.Effect <- X$Comb.Effect
Combination.name <- X$Combination.name
Sum.Mar.Effect <- X$Sum.Mar.Effect
TIE <- X$TIE
range.change <- X$range.change
}
}
else {
data.column <- data[, -1]
column.com <- list()
Range.com <- list()
for (j in 1:choose(length(data.column), order)) {
column.com[[j]] <- colnames(data.column)[c(combn(length(data.column),
order)[, j])]
if (measure == "TIE") {
Range.com[[j]] <- comb.prem(data = data, column = column.com[[j]],
Restriction = Restriction, order = order)$range.change
}
if (measure == "Order.Diff") {
Range.com[[j]] <- comb.prem(data = data, column = column.com[[j]],
Restriction = Restriction, order = order)$max.order.diff
}
}
column.com <- as.data.frame(matrix(unlist(column.com),
ncol = order, byrow = TRUE))
Compare <- as.data.frame(cbind(column.com, unlist(Range.com)))
Compare <- Compare[order(Compare[, (order + 1)], decreasing = TRUE),
]
if (measure == "TIE") {
colnames(Compare) <- c(colnames(Compare)[1:order],
"range.TIE")
}
if (measure == "Order.Diff") {
colnames(Compare) <- c(colnames(Compare)[1:order],
"max.Order.diff")
}
return(Compare)
}
if (!missing(column)) {
if (length(column) == 2) {
Treatment.coefs.name <- Treatment.Effect[, 2]
if (ncol(Treatment.Effect) >= 3) {
for (i in 1:nrow(Treatment.Effect)) {
for (j in 3:ncol(Treatment.Effect)) {
Treatment.coefs.name[i] <- paste(Treatment.coefs.name[i],
Treatment.Effect[i, j], sep = ".")
}
}
}
Treatment.coefs.name <- gsub(" ", ".", Treatment.coefs.name)
}
}
if (missing(column)) {
Treatment.coefs.name <- rownames(Treatment.Effect1)
Treatment.coefs.name <- gsub(" ", ".", Treatment.coefs.name)
}
if (missing(column)) {
Coefficients <- c()
for (i in 1:length(Treatment.coefs.name)) {
if (Treatment.coefs.name[i] %in% names(coefs)) {
Coefficients[i] <- coefs[names(coefs) == Treatment.coefs.name[i]]
}
else {
Coefficients[i] <- NA
}
}
if (outcome.type == "binary") {
Coefficients <- Coefficients/2
}
}
if (!missing(column)) {
if (length(column) == 2) {
Coefficients <- c()
for (i in 1:length(Treatment.coefs.name)) {
if (Treatment.coefs.name[i] %in% names(coefs)) {
Coefficients[i] <- coefs[names(coefs) == Treatment.coefs.name[i]]
}
else {
Coefficients[i] <- NA
}
}
if (outcome.type == "binary") {
Coefficients <- Coefficients/2
}
}
}
if (!missing(column)) {
if (length(column) == 2) {
Treatment.Effect.print <- cbind(Comb.Effect, TIE,
Combination.name)
colnames(Treatment.Effect.print) <- c("ATCE", "AMTIE",
colnames(Combination.name))
Main.matrix <- cbind(Sum.Mar.Effect, TIE, Combination.name)
colnames(Main.matrix) <- c("Sum of AMTEs", "AMTIE",
colnames(Combination.name))
Inequality.Cor <- cor(Sum.Mar.Effect, TIE)
Relative.change <- cbind(TIE, Combination.name)
colnames(Relative.change) <- c("AMTIE", colnames(Combination.name))
Relative.change <- Relative.change[order(Relative.change[,
1], decreasing = TRUE), ]
}
if (length(column) >= 3) {
Treatment.Effect.print <- cbind(Comb.Effect, TIE,
Combination.name)
colnames(Treatment.Effect.print) <- c("ATCE", "AMTIE",
colnames(Combination.name))
Main.matrix <- cbind(Sum.Mar.Effect, TIE, Combination.name)
colnames(Main.matrix) <- c("Sum of AMTEs", "AMTIE",
colnames(Combination.name))
Relative.change <- cbind(TIE, Combination.name)
colnames(Relative.change) <- c("AMTIE", colnames(Combination.name))
Relative.change <- Relative.change[order(Relative.change[,
1], decreasing = TRUE), ]
}
}
else {
Treatment.Effect.print <- Treatment.Effect1
Man.matrix <- NULL
}
Treatment.Effect.print <- as.data.frame(Treatment.Effect.print)
if (!missing(column)) {
Main.matrix <- as.data.frame(Main.matrix)
}
if (sort == TRUE) {
Treatment.Effect.print1 <- Treatment.Effect.print
Treatment.Effect.print <- Treatment.Effect.print[order(Treatment.Effect.print[,
1], decreasing = TRUE), ]
if (!missing(column)) {
Main.matrix <- Main.matrix[order(Main.matrix$Sum,
decreasing = TRUE), ]
Treatment.Effect.print <- cbind(Treatment.Effect.print[,
1:2], Treatment.Effect.print[, 3:ncol(Treatment.Effect.print)])
}
}
if (!missing(column)) {
if (sort) {
return(list(`Range of AMTIE` = range.change, AMTIE = Relative.change,
ATCE = Treatment.Effect.print, `Sum of AMTEs` = Main.matrix))
}
else {
return(list(`Range of AMTIE` = range.change, AMTIE = Relative.change,
ATCE = Treatment.Effect.print, `Sum of AMTEs` = Main.matrix))
}
}
else {
return(list(Treatment.Effect = Treatment.Effect.print))
}
}
compare1 <- INT(fit4, target.data=pred.dat, compare=TRUE, order=1) #Main effects
compare2 <- INT(fit4, target.data=pred.dat, compare=TRUE, order=2) #Two-way interactions
compare3 <- INT(fit4, target.data=pred.dat, compare=TRUE, order=3) #Three-way interactions
##Calculate barplot data
sorted.compare1 <- sort(compare1, decreasing=TRUE)
compare.names1 <- names(sorted.compare1)
sorted.compare2 <- compare2[order(-compare2[,3]),]
compare.names2 <- paste(eval(sorted.compare2[,1]), eval(sorted.compare2[,2]), sep = ":")
sorted.compare3 <- compare3[order(-compare3[,4]),]
compare.names3 <- paste(eval(sorted.compare3[,1]), eval(sorted.compare3[,2]), eval(sorted.compare3[,3]), sep = ":")
compare.names1r <- compare.names1
compare.names1r <- sub("HR", "Human Rights",compare.names1r)
compare.names1r <- sub("Rep", "Reputation",compare.names1r)
compare.names1r <- sub("Pow", "Power",compare.names1r)
compare.names1r <- sub("Gov", "Government",compare.names1r)
compare.names1r <- sub("Reg", "Region",compare.names1r)
compare.names1r <- sub("All", "Alliance", compare.names1r)
compare.names2r <- compare.names2
compare.names2r <- sub("HR", "Human Rights", compare.names2r)
compare.names2r <- sub("Rep","Reputation", compare.names2r)
compare.names2r <- sub("Pow", "Power", compare.names2r)
compare.names2r <- sub("Gov", "Govenment", compare.names2r)
compare.names2r <- sub("Reg", "Region", compare.names2r)
compare.names2r <- sub("All", "Alliance", compare.names2r)
compare.names3r <- compare.names3
compare.names3r <- sub("HR", "Human Rights", compare.names3r)
compare.names3r <- sub("Rep","Reputation", compare.names3r)
compare.names3r <- sub("Pow", "Power", compare.names3r)
compare.names3r <- sub("Gov", "Government", compare.names3r)
compare.names3r <- sub("Reg", "Region", compare.names3r)
compare.names3r <- sub("All", "Alliance", compare.names3r)
conjoint$Gov <- factor(conjoint$Gov, ordered = TRUE,
levels=c("Democracy", "Anocracy", "Autocracy"))
conjoint$Pow <- factor(conjoint$Pow, ordered = TRUE,
levels =c("top 15 most powerful", "top 25 most powerful", "top 50 most powerful"))
conjoint$Reg<- factor(conjoint$Reg, ordered = TRUE,
levels = c("Middle East North Africa", "subSaharan Africa",
"Asia", "South America", "Europe"))
conjoint$All <- factor(conjoint$All, ordered = TRUE,
levels = c("Defensive Alliance", "Offensive Alliance", "Neutrality Pact"))
conjoint$Rep <- factor(conjoint$Rep, ordered = TRUE,
levels = c("Consistently fulfills alliance commitments",
"Inconsistently fulfills alliance commitments",
"Rarely fulfills alliance commitments"))
conjoint$HR <- factor(conjoint$HR, ordered = TRUE,
levels = c("Excellent human rights record",
"Mediocre human rights record",
"Routinely violates human rights"))
conjoint$select <- factor(conjoint$selected, ordered = FALSE,
levels = c("0", "1"))
fit10 <- CausalANOVA(formula = select ~ Rep + Pow +
Reg + All + HR + Gov,
data = conjoint, diff=TRUE, pair.id = conjoint$case_no,
cluster = conjoint$random_id, nway=2, select.prob = TRUE,
seed = 45869, boot = 5000)
